Optimal multi-electrode transcutaneous stimulation with high locality and intensity

ABSTRACT

Methods, apparatus, and systems are disclosed for optimization techniques and a realistic 3D model to design optimal parameters for transcutaneous stimulation to achieve focalized stimulation of a target tissue such as the spinal cord, brain or other internal organ. The methods, apparatus, and systems include generation of a 3D model from a CT/MRI image, as well as an optimization algorithm that enables stimulation of any target location (e.g., on the dorsal root, or on the dorsal column) with any orientation at high precision.

CROSS-REFERENCE TO RELATED APPLICATIONS

This application claims priority to, and is a 35 U.S.C. § 111(a)continuation of, PCT international application number PCT/US2017/064970filed on Dec. 6, 2017, incorporated herein by reference in its entirety,which claims priority to, and the benefit of, U.S. provisional patentapplication Ser. No. 62/532,477 filed on Jul. 14, 2017, incorporatedherein by reference in its entirety, and which claims priority to, andthe benefit of, U.S. provisional application Ser. No. 62/430,490 filedon Dec. 6, 2016, incorporated herein by reference in its entirety.Priority is claimed to each of the foregoing applications.

The above-referenced PCT international application was published as PCTInternational Publication No. WO 2018/106843 on Jun. 14, 2018, whichpublication is incorporated herein by reference in its entirety.

STATEMENT REGARDING FEDERALLY SPONSORED RESEARCH OR DEVELOPMENT

Not Applicable

NOTICE OF MATERIAL SUBJECT TO COPYRIGHT PROTECTION

A portion of the material in this patent document may be subject tocopyright protection under the copyright laws of the United States andof other countries. The owner of the copyright rights has no objectionto the facsimile reproduction by anyone of the patent document or thepatent disclosure, as it appears in the United States Patent andTrademark Office publicly available file or records, but otherwisereserves all copyright rights whatsoever. The copyright owner does nothereby waive any of its rights to have this patent document maintainedin secrecy, including without limitation its rights pursuant to 37C.F.R. § 1.14.

BACKGROUND 1. Technical Field

The technology of this disclosure pertains generally to systems fortherapeutic stimulation and imaging, and more particularly tomulti-electrode transcutaneous stimulation for therapy and medicalimaging.

2. Background Discussion

Epidural spinal cord stimulation (eSCS) has been used for restoration ofmotor functions in spinal cord injury (SCI), pain management, andspasticity control. However, it is invasive that needs surgery toimplant electrodes into the body.

As an alternative, the transcutaneous spinal cord stimulation (tSCS) wasused to achieve similar effect as eSCS in a noninvasive way. It has beenshown that tSCS is able to elicit the locomotor-like movements inhealthy subjects as well as in spinal cord injury (SCI) subjects. tSCShas also been demonstrated to control pain and suppress spasticity.Conventional tSCS uses one or two large electrodes for stimulation,resulting in non-focal current flows in the spinal cord. Due to largeactivated area, it is difficult to understand the underlying mechanism,for example which part of the activated area contributes to thetreatment. Further, the two-electrode montage is not able to targetmultiple sites simultaneously, which limits its effectiveness. Whenstimulating different targets, it is often desirable to frequentlychange the electrode locations. Therefore, it becomes inconvenient forstimulating a neural network dynamically.

More recently, multiple electrode arrays comprising 3*3 electrodes, 3*7and 3*8 electrodes have been adopted for transcutaneous spinal cordstimulation, which allow multisite stimulation. Compared to single sitestimulation, the multisite stimulation has been demonstrated to inducemore effective stepping movements and higher amplitude of EMG activityin healthy subjects. However, because their stimulation parameters arejust chosen by experience, the induced current is still not focused. Inaddition, certain sensitive regions such as the bladder can not beavoided. All of these together limit their capability of modulating theneuronal circuits precisely.

To improve the focal accuracy of the stimulation, ring configuration iswidely used to enhance the focality of the epidural spinal cordstimulation (eSCS) and transcranial current stimulation (tCS). It isfeatured by an anode (cathode) electrode surrounding by four cathode(anode) electrodes. Generally, it is good at stimulating the radialorientation, and has difficulty in dealing with tangential orientation.Another limitation lies in that when the target is not underneath anyelectrode, the ring configuration is not able to stimulate the targetprecisely. In addition, it is not able to avoid certain regions either.

A more effective and focal stimulation can be achieved by preciselyconstructing a spinal cord/head model and taking advantage ofoptimization methods. Several optimization methods have been developedfor transcranial direction/alternating current stimulation (tDCS/tACS)as well as deep brain stimulation (DBS). However, there has not been anywork that uses optimization methods to obtain higher intensity or focalaccuracy for transcutaneous spinal cord stimulation. In tDCS/tACS, theconventional optimization methods either maximize intensity at thetarget, which results in very low focal accuracy of stimulation (e.g.,maximum intensity method), or maximize the focal accuracy at the expenseof low intensity (e.g., Linear Constrained Minimum Variance (LCMV)). Inaddition, in the LCMV method and it variants, a hard constraint isenforced to meet the specified intensity and orientation at the target,which may lead to infeasible solution when the specified intensity ishigh or the target region is large. In addition, the hard constraintlimits the degree of freedom of the problem, which hinders it fromobtaining a better solution with higher intensity or focality. Takentogether, an optimization method that that is able to always provide afeasible solution, as well as to optimize both intensity and focality atthe same time is highly desirable.

Another limitation of conventional optimization methods is that theyrequire the clinician to specify location and intensity at the target,which is usually unknown in most applications. Recently, a method basedon the reciprocity principle was proposed, which enables the EEG signalto be used as a guide for designing stimulation patterns withoutspecifying the location of target. However, the stimulation parameterschosen by this method is just empirical, and it is only able to dealwith simple situations, such as cases in which a single focal sourcepredominates. It is unable to handle complicated cases, such as multipletargets, spatially extended targets with different orientations indifferent parts, or containing brain regions to be avoided. In themultiple targets case, it is possible that this method may stimulate theaveraged location of these targets. In sum, a better method for guidingthe stimulation is needed.

From the hardware perspective, it is a big challenge to control a largenumber of stimulation channels independently. In addition, the key to ahighly spatially focused stimulation is enabling precise control ofcurrents in each stimulation channel. Enabling necessarily high count ofindependently controlled current sources is challenging.

Thus, in order to use an optimization algorithm to its highestprecision, it is highly desirable to have a high channel-countstimulator with precise output current parameters for each channel andseemingly immediate dynamic updates to those parameters in real-time.

BRIEF SUMMARY

One aspect of the present description is a system and method forgeneration of a 3D model from a CT/MRI image, as well as an optimizationalgorithm that enables stimulation of any target location (e.g., on thedorsal root, or on the dorsal column) with any orientation at highprecision. In addition, the systems and methods of the presentdescription are capable of stimulating not only single target but alsomultiple targets, as well as avoid certain regions within the anatomy ofinterest. For hardware, the system provides a hardware platform withcapabilities to enable the optimization algorithm to its highestprecision: a high channel-count stimulator with precise output currentparameters for each channel and seemingly immediate dynamic updates tothose parameters in real-time.

One aspect of the present description is a method for optimal and focaltranscutaneous spinal cord stimulation.

Another aspect of the present description is a method for transcranialcurrent stimulation.

The systems and methods disclosed herein provide better results thanother state-of-the-art methods in terms of directional intensity andfocality, and can be extended to stimulation of other internal organs,including the brain.

Further aspects of the technology described herein will be brought outin the following portions of the specification, wherein the detaileddescription is for the purpose of fully disclosing preferred embodimentsof the technology without placing limitations thereon.

BRIEF DESCRIPTION OF THE SEVERAL VIEWS OF THE DRAWING(S)

The technology described herein will be more fully understood byreference to the following drawings which are for illustrative purposesonly:

FIG. 1 shows a perspective view of a 3D spinal cord model generated froman individual MRI/CT image in accordance with the present description.

FIG. 2 shows a perspective view of the spinal cord model of FIG. 1 .integrated with a multi-electrode stimulation array.

FIG. 3 shows a perspective view of the spinal cord model of FIG. 1 aftermeshing.

FIG. 4A shows an image of a spinal column with a single simulatedtarget.

FIG. 4B shows an image of a spinal column with a multiple simulatedtargets.

FIG. 4C shows an image of a spinal column with a single simulated targetone nearby avoidance region.

FIG. 5A through FIG. 5C show the stimulation results for theoptimization method of the present description, with FIG. 5A showing theelectrode weighting (stimulation parameter at each electrode), FIG. 5Bshowing the intensity E-field at the white matter, and FIG. 5C showingthe Directional E-Field in the desired direction (y-axis).

FIG. 6A through FIG. 6C show the stimulation results for theoptimization method of the present description for target orientationalong the z-axis (tangential to the electrode), with FIG. 6A showing theelectrode weighting (stimulation parameter at each electrode), FIG. 6Bshowing the intensity E-field at the white matter, and FIG. 5C showingthe Directional E-Field intensity in the desired direction (z-axis).

FIG. 7 shows an image of a 3D head model including scalp, skull, CSF,and cortex.

FIG. 8 shows a 3D head model with electrodes and course mesh pattern.

FIG. 9 shows a 3D head model after meshing is applied.

FIG. 10A shows a top view of a 3D brain model with mesh pattern and thetarget locations.

FIG. 10B shows a side view of the 3D brain model with avoidance region.

FIG. 11A shows the resulting electrical field for the LCMV method.

FIG. 11B shows the resulting electrical field for the optimizationmethod of the present description.

FIG. 12A shows the weighting for each electrode used for the LCMVmethod.

FIG. 12B shows the weighting for each electrode used for theoptimization method of the present description.

FIG. 13 shows a plot of the focality of the LCMV method, maximumintensity method and the optimization method of the present descriptionwith varying electrode number.

FIG. 14 shows a schematic diagram of a stimulation system withmulti-channel stimulation SoC in accordance with the presentdescription.

FIG. 15 shows an image of E-field from in-vitro stimulation usingsimultaneous multi-channel stimulation.

FIG. 16 shows a diagram of an electrode array using five electrodes formulti-channel stimulation.

FIG. 17 shows a flow diagram for construction of a spinal cord 3D model.

FIG. 18 shows a schematic diagram of a transcutaneous multi-electrodearray and stimulator in accordance with the present description.

FIG. 19 shows a schematic diagram of a therapeutic spinal cordstimulation system positioned on a patient.

FIG. 20 is a schematic block diagram of a wireless sensing andstimulation system with a controller configured to be held by thepatient or the clinician.

FIG. 21 is a schematic block diagram of a sensor network in accordancewith the present description.

FIG. 22 shows a diagram for an exemplary modulator and modulator outputfor the system and methods of the present description.

DETAILED DESCRIPTION

The key challenges for noninvasive transcutaneous spinal cordstimulation and transcranial current stimulation lie in the ability toprovide high spatial and temporal resolution, with a high degree offocal accuracy, while using the correct intensity and directionalityusing external current injection though an electrode on the skin.

To accommodate the above mentioned challenges, a multi-electrode systemis described below. The multi-electrode system of the present technologyincludes novel hardware, software, algorithms, and user interface thatenable an experiment methodology and development of clinical protocolsby (1) a computational 3D model based on MRI/CT data to facilitate thetreatment, efficacy, and safety of the device; (2) a novel sourcelocalization algorithm to obtain precise target localization and focalprecision; (3) optimal stimulation patterns with desired intensityguided by a novel source localization algorithm; (4) miniaturized andfully integrated hardware electronics to independently drive each of themulti-electrodes, and thus offering the portability and flexibility atof a lower cost system; and (5) a sophisticated but user-friendlygraphical user interface that will enable home use.

The system and methods of the present description are able to achievefocal and precise stimulation for transcutaneous spinal cordstimulation, and can be used to target single or multiple sites, as wellas avoid certain regions. It is also capable of targeting any locationwith any orientation depending on different applications. Optimizationtechniques are combined with precise spinal cord modeling to provideoptimal stimulation parameters for transcutaneous spinal cordstimulation. The system and methods of the present description makesfull use of each electrode in a multi-electrode array to achieve focaland precise stimulation. It has also been demonstrated that themethodology of the system and methods of the present description can beequivalently applicable to transcranial current stimulation, andpossibly the noninvasive internal organ stimulation.

To overcome the limitations of existing systems, the system and methodsof the present description employ a novel optimization technique thatprovides a solution with both high intensity and high focal accuracywithin the safety constraints. In particular, for the spatially extendedtarget, the algorithm of the system and methods of the presentdescription is able to provide various current intensity distribution(e.g., uniform, smooth, Gaussian, etc.) at the target depending on theapplication.

Another limitation of previously disclosed optimization methods is thatthey require the clinician to specify location and intensity of thetarget, which is usually unknown in most applications. With use ofdynamic EEG brain imaging disclosed in PCT International Application No.PCT/US2016/050452 filed on Sep. 6, 2016 and published as WO 2017/044433A1 on Mar. 16, 2017, and precise EEG source localization, the systemsand methods of the present description is able to provide the accurateinformation of the target location, number as well as orientation, so asto enable a precise stimulation.

The EEG brain imaging system of the present description is able toprovide much higher temporal resolution in the range of millisecondsrather than seconds. In addition to guide dynamic stimulation of theneural networks, the concurrent EEG brain imaging will also offerreal-time feedback of the neuromodulation. Thus, a closed-loopstimulation is contemplated. Similarly, for transcutaneous spinal cordstimulation, EMG inverse imaging reconstructed from the electricalpotential recorded by the surface electrodes may be used as a guidance.

1. Mathematical Formulation and Modeling of Spinal Cord Model

FIG. 1 shows an image of a computer generated spine model 10 generatedfrom a CT/MRI image 12. A spinal cord model with different tissues isassumed, with each tissue having isotropic conductivity. Assuming Nelectrodes on the skin of the back, and one large return electrode onthe belly. We use a vector x ∈

to represent the injected current at each electrode. Further, wediscretize the whole spinal cord volume into M voxels, and use vector e∈

to represents the electrical field of each voxel resulting from theelectrical stimulation. Since we further consider the orientation of theelectrical field, so the vector e has a dimension of 3M*1. Underquasi-static condition, the electrical field e in the voxels and thestimulation parameters x at the electrodes has a linear relationship:e=Kx.  Eq. 1

Here the coefficient matrix K is called “lead field matrix”, whichdescribes the one to one mapping between each electrode and each voxel.Specifically, the (i,j)^(th) entry of K denotes the electrical field atthe i^(th) voxel due to an unit current stimulation at the j^(th)electrode. Therefore, the electrical field in each voxel is a linearsuperposition of that resulting from the injected current from eachelectrode. The K matrix can be calculated by constructing a spinal cordmodel and solving the Maxwell's equations with the boundary elementmethod (BEM) or finite element method (FEM).

To calculate the lead field matrix, a spinal cord model is firstconstructed. A 3D spinal cord model is generated based on ahigh-resolution CT/MRI image, which includes the following steps: imagesegmentation, electrode model construction, meshing. The FEM method tocalculate the lead field matrix in Eq. 1.

(a) Image Segmentation

Referring to FIG. 1 , the torso 14 of model 10 is segmented intodifferent tissues according to the gray level on the MRI/CT image 12.The segmentation can be done manually in software such as Solidworks, orautomatically in software such as MeVisLab. For tissues that aredifficult to be identified from the image 12, portions of the model maybe built manually using software such as Solidworks, COMSOL, or thelike.

After image segmentation, 2D segmentation results are converted into a3D tissue model for each tissue, e.g., using Solidworks or likesoftware. All of the tissue are then assembled together to form a whole3D model 10 as shown in FIG. 1 . FIG. 1 shows a 3D model 10 withdifferent tissues, such as skin 16 (including stratum corneum (SC),stratum germinativum (SG) and dermis layers), fat/muscle 18, vertebrae20, spinal cord/nerve 22 (gray and white matters), cerebrospinal fluid(CSF), etc.

(b) Electrode Model Construction

Referring to FIG. 2 , after constructing the 3D tissue model 10 for thespinal cord, we import the model into software such as COMSOL, andconstruct an electrode model 30. The electrode model 30 includes amulti-electrode array of electrodes 32 on the back and one returnelectrode 33 on the belly. FIG. 2 shows an array 30 of 9 (rows)*7(columns) of small electrodes 32 for stimulation and one large electrodeas the return 33. In one exemplary model, the electrodes 32 are composedof a stainless steel material with thickness of 60 μm and diameter 1 cm.It is appreciated that other electrode array configurations, as well aselectrode composition/dimension, may be contemplated. The electrodematerial, size, pitch, thickness, position, etc. may be modified asappropriate in the modeling software.

(c) Meshing

Referring to FIG. 3 , the whole model (which includes the tissue model10 and electrode model 30) is then discretized into a large number ofvoxels to form a finite element model 34.

(d) Lead Field Matrix Calculation

To calculate the lead field matrix of Eq. 1, any finite element methodsoftware, such as COMSOL, may be implemented. For each electrode 32, aunit current (density) assigned, and the resultingvoltage/E-field/current (density)/activation function is calculated ateach voxel. The obtained values form a vector, which will become thecorresponding column for the lead field matrix. Repeating this processwill result in the whole lead field matrix. Table 1 lists theconductivity values for different tissues used to calculate the leadfield matrix.

2. Optimization Methods

This section details a novel method to design and configure optimalstimulation parameters for transcutaneous spinal cord stimulation (tSCS)to improve the focal accuracy of tSCS. This method overcomes thelimitation of maximum intensity and LCMV methods, and is able to obtainboth high intensity and focality at the target.

(a) Safety Limit

First, the safety limit of transcutaneous spinal cord stimulation isdetermined. With a carrier frequency of 10 kHz, the subject can easilytolerate a stimulation intensity from 30 mA to 200 mA and not feel pain.Therefore, a reasonable safety criteria assumes a current intensity foreach electrode that is limited to 100 mA, with the total currentintensity limited to 200 mA, and the intensity in the avoidance regionis at least 10 times smaller than that in the target region. UsingI_(max) to represent the maximum current at each electrode, I_(total) todenote the maximum total current injected to the body, and ratio torepresent the intensity ratio between the target and avoidance region,the safety criteria are summarized according to FIG. 2 through FIG. 5 :

$\quad\left\{ \begin{matrix}{{{x_{i}} \leq I_{{ma}\; x}},\;{i = 1},\ldots\;,\; N} & {{{Eq}.\mspace{11mu} 2}\mspace{14mu}} \\{{\sum\limits_{i = 1}^{N}{x_{i}}} \leq {2*I_{total}}} & {{{Eq}.\mspace{11mu} 3}\mspace{14mu}} \\{{\sum\limits_{i = 1}^{N}x_{i}} = 0} & {{{Eq}.\mspace{11mu} 4}\mspace{14mu}} \\{{Intensity}_{avoid} \leq {\frac{1}{ratio}{Intensity}_{target}}} & {{{Eq}.\mspace{11mu} 5}\mspace{14mu}}\end{matrix} \right.$

Note that since the total positive current is equal to the totalnegative current (Eq. 4), therefore the sum of absolute value of thecurrent should be twice of the total current injected to the body (Eq.3). It is worth noting that the suggested safety limit valuesI_(max)=100 mA and I_(total)=200 mA can be easily modified for differentapplications.

(b) Optimization Model

The key challenges for noninvasive stimulation technologies lie in thecapability of providing precise stimulation with both high focalaccuracy and intensity in the desired direction. Current optimizationmethods either maximize the intensity at the target by satisfying thefocal accuracy (e.g., maximum intensity method), or maximize the focalaccuracy at the expense of low intensity (e.g., Linear ConstrainedMinimum Variance (LCMV)). In addition, in the LCMV method and itvariants, a hard constraint is enforced to meet the specified intensityand orientation at the target, which may lead to infeasible solutionwhen the specified intensity is high or the target region is large. Toovercome these limitations, the optimization method of the presentdescription always provides a feasible solution that includes both highintensity and focal accuracy within the safety constraints. This modelcan be expressed as follows:

$\begin{matrix}{{x_{proposed} = {{\arg\mspace{11mu}{\min\limits_{x}{\frac{1}{w}{{Kx}}_{2}^{2}}}} - {\lambda*e_{0}^{T}{Cx}}}}{{subject}\mspace{14mu}{to}\mspace{14mu}\left\{ \begin{matrix}{{{x_{i}} \leq I_{\max}},{i = 1},\ldots\;,N} \\{{\sum\limits_{i = 1}^{N}{x_{i}}} \leq {2*I_{total}}} \\{{\sum\limits_{i = 1}^{N}x_{i}} = 0} \\{{Intensity}_{avoid} \leq {\frac{1}{ratio}{Intensity}_{target}}}\end{matrix} \right.}} & {{Eq}.\mspace{14mu} 6}\end{matrix}$where the constant w is equal to the ratio between the total number ofvoxels and the number of targeted voxels. The first term is the focalityterm, the second term is the intensity on the desired direction. Theparameter λ balances these two objectives and controls the relativeimportance of the focality and directional intensity, e.g., by applyinga weighting to the focality and directional intensity terms. It can beset empirically or automatically by the L-curve method or crossvalidation method. This optimization problem is convex, which can besolved by software such as CVX efficiently.

The algorithm detailed in Eq. 6 is able to deal with any target locationon any tissue (e.g., bone, white matter) with any target orientation. Itcan not only deal with multiple targets, but also avoid activatingcertain sensitive regions. In addition, it is very easy to incorporatevarious safety constraints into the model. Rather than setting a hardconstraint at the target like the LCMV method, it allows a range of theintensity by changing the parameter λ. It optimizes intensity andfocality simultaneously.

To further improve the focality, we can replace the L2 norm with L1 normto impose sparsity on the stimulated area, as follows

$\begin{matrix}{x_{proposed} = {{\arg\mspace{11mu}{\min\limits_{x}{\frac{1}{w}{{Kx}}_{1}}}} - {\lambda*e_{0}^{T}{Cx}}}} & {{Eq}.\mspace{14mu} 7}\end{matrix}$

Furthermore, for the spatially extended target, the algorithm is able toprovide various current intensity distribution (e.g., uniform, smooth,Gaussian, etc.) at the target depending on the applications. This can beachieved by imposing sparsity on a transform domain, rather than theoriginal domain:

$\begin{matrix}{x_{proposed} = {{\arg\mspace{11mu}{\min\limits_{x}{\frac{1}{w}{{D({Kx})}}_{1}}}} - {\lambda*e_{0}^{T}{Cx}}}} & {{Eq}.\mspace{14mu} 8}\end{matrix}$where the operator D is the total variation operator (first orderdifference operator) if uniform distribution is desired; or Laplacianoperator if smooth distribution is desired. If highly smooth/Gaussiandistribution is desired, total generalized variation (TGV) may be usedfor the second term.

(c) EEG/EMG Inverse Image Guided Optimal Stimulation

One of the limitations of conventional optimization methods is that theyrequire the clinician to specify the location and intensity of thetarget area, which is unknown in most applications.

Different from the reciprocity principle-based methods that use an EEGsignal as a guide, the system of the present description may beconfigured to use a dynamic EEG brain image system (as described inUC16-151-2FP, incorporated herein by reference in its entirety) as aguide for electrical stimulation. Compared to fMRI, EEG brain imagingprovides much higher temporal resolution (˜ms). It is able to providenot only the number, the location but also the orientation of targets.With precise EEG source localization, the system is able to provide theability to deal with complicated target/avoidance configurations, andwith spatial high accuracy. The concurrent EEG brain imaging not onlyprovides a guide for dynamic stimulation of complicated neural networks,but also offers real-time feedback of the neuromodulation. Thus, aclosed-loop stimulation is eminently possible.

Similarly, for transcutaneous spinal cord stimulation, the surfacemulti-electrode array may be used to record the electrical signal, andcalculate the solution of the inverse problem to find out the targets,then use the targets as a guidance for the optimization solver.

(d) Simulation Protocol

Referring to FIG. 4A to FIG. 4C, a study was performed by simulatingthree different target configurations to test the performance of ourmethod in different cases, including 1) single target case (FIG. 4A), inwhich there is only one target region 42; 2) multiple targets case (FIG.4B), in which there are two target regions 42 a and 42 b; and 3)avoidance case (FIG. 4C), in which one target 42 is stimulated alongwith a nearby avoidance region 44. For the transcutaneous spinal cordstimulation, the default target tissue is white matter. This is becausein many applications, such as locomotor behavior regulation and paincontrol, the target region is dorsal root or dorsal column. Except forwhite matter, the case of the target location on the vertebra was alsotested (e.g., a spot on the vertebra where there is a cavity was chosenso that the current can flow in to spinal cord). The results werecompared with that of white matter. For the target orientation, thedefault orientation is along y-axis (radial to the electrode), but theorientation of z-axis was also tested, and a comparison between them wasmade.

(e) Quantitative Metrics

In order to evaluate the performance of various methods quantitatively,we use the following criteria:

(i) Intensity (unit: V/m), defined as the average intensity at thetarget voxel;

(ii) Target Error (TE, unit: cm or mm), defined as the distance betweencenter of mass of the activation area and the target; and

(iii) Focality (unit: cm or mm), defined as the radius within which theaccumulative energy (square of intensity) is half of the total energy.

Considering that the target orientation also plays an important role onthe stimulation effectiveness, the following directional criteria arefurther adopted:

(i) Directional Intensity (DI, unit: V/m), defined as the averageintensity in the desired direction;

(ii) Directional target error (DTE), which is similar to target error(TE). except for using directional intensity instead of total intensitywhen calculating TE; and

(iii) Directional focality (DF, unit: cm or mm), which is similar tofocality except for using directional intensity rather than totalintensity.

In one embodiment, spinal cord stimulation parameters may comprise acurrent intensity from each electrode ranges from 0.01 mA to 250 mA,pulse width ranges from tens of μs tens of ms, and the time separationbetween stimulus or a group of stimuli ranges from 10 s to 0.001 s. Theoverall transient current intensity from all electrode is less than 250mA.

3. Spin Model Results

In this section, the qualitative and quantitative results of variousmethods for transcutaneous spinal cord stimulation are shown. In orderto make a fair comparison among different methods, the same safetycriteria (Eq. 2-Eq. 5) was applied for all of the methods.

The method of the present description was compared with several existingmethods in the literature: single large electrode, single smallelectrode, maximum intensity, ring configuration, and LCMV. FIG. 5Athrough FIG. 5C shows the stimulation results for the method of thepresent description, with FIG. 5A showing the electrode weighting, FIG.5B showing the intensity E-field, and FIG. 5C showing the DirectionalE-Field.

For single large electrode, single small electrode, and the maximumintensity methods, the intensity of E-field at the target is very high.However, the E-field is spread out, resulting in very low focalaccuracy. Among all the methods, the maximum intensity method providesthe highest directional intensity and the worst focality. In comparison,the results of ring configuration, LCMV and optimization method of thepresent description are much more focal. Compared to the other twomethods, the optimization method of the present description provides thebest results in terms of both intensity and focality.

The optimization method of the present description is able to target anylocation including any tissue such as bone and white matter. Deepersource is generally more difficult to target, in terms of intensity andfocality. For example, the results of targeting bone vs. targeting whitematter were also compared. Results show that targeting bone providesmuch higher target intensity and focality accuracy than targeting whitematter, which is reasonable since when targeting the white matter, theelectrical field is weakened by the high-resistivity bone.

Except for target locations, optimization method of the presentdescription is also able to target any orientation. FIG. 5A through FIG.5C show results we set the target orientation to be along the y-axis(radial to the electrode). FIG. 6A through FIG. 6C, the results oftarget orientation along the z-axis (tangential to the electrode) areshown, and compare the performance of different stimulation methods.FIG. 6A shows the stimulation parameters for each method, and one cansee that the weighting for each electrode is very different from that ofy-axis orientation (FIG. 5A).

For target with z-axis orientation, generally the stimulation pattern isto place an anode on one side of z-axis and a cathode on the other side.The ring configuration has difficulty in dealing with tangentialtargets, which resulting in low target intensity and focality. Toimprove the focal accuracy, the LCMV method and the optimization methodof the present description place multiple anodes with differentweightings rather than only using single anode/cathode. Compared toother methods, the optimization method of the present descriptionprovides the best focal accuracy.

The optimization method of the present description is also able to dealwith both single target and multiple targets. The optimization method ofthe present description provides better results than LCMV in terms oftarget intensity, localization accuracy, and focal accuracy. Compared tothe results of single target, both target intensity and focality arelower in the multiple targets case, due to additional constraints on thesolution.

With respect to avoidance regions, the ring pattern is unable to avoidcertain regions. In contrast, optimization methods including theoptimization method of the present description are able to achieveavoidance by constraining the intensity in the avoidance region to xtimes (e.g., x=10) lower than that in the target region. The performanceof the LCMV method and the optimization method of the presentdescription were compared in performance with dealing with an avoidanceregion close to the target region. The result shows that the LCMV avoidsthe region with the expense of shifting the activated area away from thetarget region, leading to larger target error and lower focality. Theoptimization method of the present description demonstrates much higherlocalization accuracy and focal accuracy than LCMV in this case.

The intensity and focality at the target have an inherent trade-off. Inthe optimization method of the present description, the parameter λcontrols the relative importance between the target intensity andfocality, and is critical to the stimulation results. When λ is large,more weight is put on the directional intensity term, therefore theintensity of the results will be high, and vise versa. By setting λ tobe a very small value (e.g., 0.01), we can estimate the upper bound ofthe focality. On the other hand, when setting λ to be a very large value(e.g., 1000), we can obtain an upper bound of the intensity. As λ goesto infinity, the optimization method of the present descriptionessentially becomes the maximum intensity method. To obtain a bestresult, selecting an appropriate value of λ is very important.

The LCMV method, the maximum intensity method and the optimizationmethod of the present description were compared with different λ. Theresults show that as λ becomes larger, the directional focalitydecreases while the directional intensity increases. The result ofλ=0.01 estimates the upper bound of directional focality, which is 0.59cm, and that of λ=1000 estimates that of directional intensity, which is27.3V/m. One can see that the results of λ=1000 matches with that ofmaximum intensity method. In addition, our method always obtains betterresults than LCMV in terms of focality and intensity. For example, whenλ=2, our method provides a similar intensity to LCMV, but its focalaccuracy is higher; when λ=4, it provides similar focality to LCMV, butthe intensity is much higher.

By changing λ, we can estimate an upper bound of focality and intensity.When λ is very large, we get the best intensity; when λ is very small,we get the best focality. With an appropriate parameter λ, we can obtainan elegant solution with both high intensity and focality. The parameterλ can be selected manually by experience, or automatically using L-curveor cross validation methods.

Note that the above systems and methods may be applied not only tohumans, but also to animals, e.g., rat, monkey, pig. In addition, it canbe applied to any part of the spinal cord, including cervical, thoracic,lumbar, and sacral.

4. Brain Stimulation

In this section, the systems and methods detailed above were applied totranscranial current stimulation (tCS). We demonstrate that our methodalso outperforms other state-of-the-art methods in tCS. The proceduresfor extending to noninvasive internal organ stimulation are alsodescribed briefly.

(a) 3D Head Model

Using a similar method as the spinal cord model construction, a 3D headmodel and source model can be constructed from a high-resolution MRimage. FIG. 7 shows a 3D head model built with various tissues: scalp52, skull 54, CSF 56, and cortex 58.

In addition, an electrode model containing multiple electrodes wasconstructed, where the electrode locations were devised from theinternational standard 10-10 system. FIG. 8 shows a 3D head model 60with 64 electrodes and course mesh pattern 62, and FIG. 9 shows a 3Dhead model 66 after meshing is applied. In the electrode models of FIG.8 and FIG. 9 , the electrodes use a gel material and 2 mm thick copper,with a diameter of 1.2 cm. Note that the electrode array andconfiguration may be varied by number, material, size, pitch, thicknessand position can be easily modified.

With the generated 3D head model, FEM is used to calculate the leadfield matrix in COMSOL. For the conductivity of each tissue, publishedvalues were used.

(b) Safety Limits for Transcranial Current Stimulation

For transcranial direct current stimulation (tDCS), there have beendifferent proposed safety criteria, including limits for total current,current density, charge density, duration, etc. A common criteria sharedby most literature is that the total current should be less than 2 mA.In the present study, both total current and current were restricted forindividual electrodes. Specifically, the following constraints wereused:

$\quad\left\{ \begin{matrix}{{{x_{i}} \leq I_{\max}},{i = 1},\ldots\;,N} \\{{\sum\limits_{i = 1}^{N}{x_{i}}} \leq {2*I_{total}}} \\{{\sum\limits_{i = 1}^{N}x_{i}} = 0} \\{{Intensity}_{avoid} \leq {\frac{1}{ratio}{Intensity}_{target}}}\end{matrix} \right.$

where the I_(max) is set to be 1 mA, the I_(total) is set to be 2 mA,and the intensity ratio between the target and avoidance region is setto be 10 in this study. Note that these values can be easily modified tosatisfy the requirements of different applications.

To test the performance of the optimization method of the presentdescription on transcranial current stimulation, several different caseswere simulated, including: single targets, multiple targets andavoidance region on the cortex surface. FIG. 8 shows a 3D head model 60with 64 electrodes and course mesh pattern 62, and FIG. 9 shows a 3Dhead model 66 after meshing is applied. FIG. 10A shows a top view of a3D brain model 70 with mesh pattern and the target locations 74. FIG.10B shows a side view of the 3D brain model 70 with avoidance region 76.For multiple targets 74, two symmetric points on the left and righthemispheres were used. For the avoidance region 76, part of the auditorycortex was used.

5. Brain Model Results

The performance of the optimization method of the present descriptionwas tested and compared to various state-of-the-art methods, including:single small electrode, ring pattern, maximum intensity, weighted leastsquare (WLS), and LCMV. The same safety limits (including total currentand current for individual electrode) were applied to the threeoptimization methods: WLS, LCMV and the optimization method of thepresent description. For the other three methods that do not useoptimization, only the total current was restricted. Additionally, theLCMV method and optimization method of the present description were alsocompared in the following cases: single target, multiple targets, andsingle target with an avoidance region.

FIG. 11A shows the resulting electrical field for the LCMV method andFIG. 11B shows the resulting electrical field for the optimizationmethod of the present description. FIG. 12A shows the weighting for eachelectrode used for the LCMV method and FIG. 12B shows the weighting foreach electrode used for the optimization method of the presentdescription. The optimization method of the present description showsbetter intensity and focality than LCMV (and WLS—not shown). The singleelectrode and maximum intensity methods (not shown) had high intensityat the target, but the focality of the electrical field was very low.The ring configuration (not shown) had a relatively focal activationarea around the target, but weak intensity compared to that on thetemporal cortex.

For one single prominent target, the optimization method of the presentdescription performed better than LCMV in terms of intensity,directional intensity (DI), target error (TE) and focality.

For multiple targets, the two symmetric points 74 on the left and righthemispheres (FIG. 10A) were set to be the targets. While both methodswere able to target both sites, and the focal accuracy is lower thanthat in the single target case. Compared to LCMV, the optimizationmethod of the present description shows higher focal accuracy.

As illustrated in FIG. 10B, part of the auditory cortex is set to be theavoidance region 76. To avoid the auditory cortex, the activated area ofLCMV shifted to the top, resulting in larger target error and lowerfocal accuracy. In comparison, the optimization method of the presentdescription provides a more focal result with higher localizationaccuracy.

The influence of the electrode number on the focal accuracy of thestimulation was also studied. LCMV, maximum intensity and theoptimization method of the present description were compared withelectrode numbers of 21, 32 and 64. FIG. 13 shows a plot of the focalityof LCMV method, maximum intensity method and the optimization method ofthe present description with different electrode number. One can seethat as the electrode number increases, the focality improves. Inparticular, the optimization method of the present description providesthe best focal accuracy.

In sum, the results showed that the optimization method of the presentdescription is able to target any location with any orientation that iseither specified by the user or provided by the EEG-brain imaging. It isable to deal with different complicated cases, such as single target,multiple targets, and target with avoidance region. It is worth notingthat the result of single target is the best, since it has lessconstraints hence has a higher degree of freedom to select a bettersolution. Compared to the maximum intensity and LCMV methods, theoptimization method of the present description maximizes the intensityand focality at the same time and achieves better focality andintensity. The benefits of the optimization method of the presentdescription are especially significant when a certain region is not tobe activated. The optimization method of the present description is ableto avoid certain region without shifting the activated area too much,while the LCMV shifts the activation area, resulting in lowerlocalization accuracy.

Results show that increasing electrode number helps to improve thestimulation. Therefore, an exemplary system may incorporate a densemulti-electrode array with up to 256 electrodes to achieve highstimulation precision. An exemplary configuration of an ultra-denseelectrode array for brain imaging may be found in PCT InternationalApplication No. PCT/US2016/050452 filed on Sep. 6, 2016 and published asWO 2017/044433 A1 on Mar. 16, 2017, incorporated herein by reference inits entirety.

6. Internal Organ Stimulation

The methodology and procedures detailed above may also be applied tostimulate an internal organ (e.g., stomach, intestine, colon, spleen, orthe like organ) noninvasively. The first step would be to construct arealistic 3D model for the internal organ (including other parts oftorso outside of the internal organ) based on high-resolution structuralimages such as MRI and CT. In addition, a multiple electrode model onthe skin should also be built in software such as but not limited toSolidWorks. Then the models can be imported into FEM software such asbut not limited to COMSOL, where the models are meshed into large numberof voxels, and the lead field matrix are calculated. After that, we canapply the same optimization methods to stimulate certain regions on theinternal organ with high localization and focal accuracy. Concurrently,we can use the stimulation electrodes for recording purpose(same-electrode stimulation and recording) to provide feedback for thestimulation.

7. Hardware System

Commercially available or existing stimulators are bulky due to theirimplementation using off-the-shelf components and wired connections. Inaddition, they lack the capability of performing simultaneousstimulation and recording to unravel the complex dynamics beneath oradjacent to the stimulating electrode. To overcome the currentlimitations, the stimulation system 100 shown in FIG. 14 uses aminiaturized, versatile and flexible wireless stimulator (e.g.,stimulator 200 shown in FIG. 17 ) in which each channel can beindependently programmed.

In addition to the optimization algorithm, the key to a highly spatiallyfocused direct current stimulation is enabling precise control ofcurrents in each stimulation channel. Thus, the optimization algorithmdetailed above may be used to its highest precision by providing ahardware platform as shown in stimulation system 100 of FIG. 14 that hasthe capability to enable this precision: a high channel-count stimulatorwith precise output current parameters for each channel and seeminglyimmediate dynamic updates to those parameters in real-time.

System 100 addresses these needs by offering an advantage of fastindependent control of stimulation parameters for each of the electrodes102 in stimulation array 104 (FIG. 14 shows 64 electrodes 102, howeverthe array 104 may ideally use 256 or more electrodes), as well assimultaneous recording of neural signals from the same electrodes with anovel algorithm which allows immunity to stimulation artifacts that willotherwise corrupt the neural signals. In addition, the system 100includes a two-way wireless data link 116, and a graphical userinterface 118 to control the system wirelessly by means of a laptop,tablet or mobile phone, which also executes the optimization algorithmsvia application programming instructions for real-time stimulationconfiguration updates. All stimulator system 100 components arepreferably integrated for ease of operation, and contain respectivesafety features allowing safe animal and human in-vivo testing.

The stimulator system 100 includes a stimulator System-On-Chip (SoC)110, for implementation of a modulation system which supports thestimulator system 10. This SoC stimulator 110 is capable of drivingstimulation current into electrode-tissue interfaces with wide ranges ofimpedance values. It is important to point out that high-densityelectrode implies smaller electrode size while small electrodeinevitably results in higher electrode-tissue interface impedance. Thus,in order to accommodate this, special high voltage CMOS devices areemployed as current sources in the design of the stimulator, resultingin the current driver's operating voltage range of +/−15V.

Another functional advantage of stimulator 110 is the capability forfast programming of each stimulation channel independently. Thestimulator 110 includes supply regulators and a data receiver 122, andtwo layers of digital controllers on board, with one global controller124 and multiple local controllers 128, both of which accepts datapackets for stimulation configuration over a single input at a fast rateof 2 Mb/s. This in turn allows the ability to digitally program everychannel with individual current settings at once within a fewmilliseconds of time. This capability is implemented to facilitatereal-time updates to stimulation based on neural signal feedback andoptimization algorithm during the stimulation session. The miniaturizedform factor of the multi-channel stimulator chip 110 allowsimplementation of the stimulator system 10 in a compact packaging,usable for both clinical or personal (at home) applications.

One or more of the FPGA114, controller 124, or external processingdevice (e.g., tablet operating GUI 118) may comprise applicationprogramming 130 stored in memory 134 for executing the generating theoptimized stimulation parameters and/or providing the stimulationcommands delivered to the electrode array 104.

Moreover, the multi-channel stimulator SoC 110 may be configured toevoke a 2-D electric field pattern in-vitro (see FIG. 16 and Table 2).This is accomplished by selectively choosing the stimulation channels(e.g., channels 1, 9, 3 and 11 in FIG. 16 ), intensities, polarities,and return channel (electrode 6 in FIG. 16 ) in an the electrode array104). FIG. 15 shows an image of recorded E-field from in-vitrostimulation using simultaneous multi-channel stimulation, as describedabove.

The on-board FPGA 114 (which may also be a microcontroller (MCU))functions as the controller between recording amplifiers 106 (withmultiplexer 108 and ADC 112), stimulator SoCs 110, and the wireless link116. Custom firmware is included for an efficient communication of datapackages between the components. The implementation of the two-waywireless link 116 and GUI software 118 leverage a neural interfaceapproach. GUI 118 may also include a safety control which allows asingle button stimulation turn off command. Finally, the full system 10integrates the stimulator 110 and recording sections with a rechargeablebattery (not shown) into a portable module with a form factor and sizecomparable to that of a smartphone. The software can be installed on alaptop or a mobile phone and use Bluetooth Low Energy (BLE) protocol toestablish the two-way link with the hardware module at a rate of 240kbps.

The technology previously described herein can be used to evaluatediseases associated with movement disorder (e.g., spinal cord injury,stroke, cerebral palsy, and Parkinson's disease) or regulation of thepatient autonomic system (e.g., blood pressure, gastrointestinalmotility, and inflammatory responses). For example, FIG. 17 shows a flowdiagram for a method 150 for construction of a 3D spinal cord model. Dueto the fact that the shape of human body varies from patient to patient,the stimulation parameters must be personalized to achieve optimaltreatment efficacy without causing tissue damage. Building a 3D spinalcord model is thus a critical component for the use of themulti-electrode array to achieve high focality and intensity stimulationwith optimized stimulation parameters using the methodology previouslydescribed herein. In the model construction, the patient first takesMRI/CT image at step 152 to derive his/her MRI/CT image. Subsequently, amechanical model based on MRI/CT images is built at step 154, preferablyincluding cerebrospinal fluid (CSF), spinal cord (gray and whitematters), vertebrae, muscle, skin (e.g., stratum germ inativum andstratum corneum), as well as fat tissues. Each physical layer would thenbe assigned with corresponding electrical and thermal property (e.g.,conductivity and permittivity) for simulation. Multi-physic simulationsoftware (e.g., COMSOL) can, for example, be used to import the modelfor EM and thermal simulation to estimate the boundaries of safestimulation parameters to avoid tissue damage. The 3D model from step154 is incorporated in the optimization of the stimulation parameters atstep 156, based on the set focality and intensity at the target or theset no-stimulation region in the spinal cord using the method previouslydescribed herein. This modeling flow can be used not only for spinalcord modeling, but also for the modeling of other organs of human body,such as brain, stomach, intestine, colon, nerves, heart, bladder, and soon. It is also not limited to transcutaneous electrical stimulation, butcan be applied to other non-invasive or invasive neuromodulation scheme,such as focused ultrasound or implantable neural stimulator.

Conventional transcutaneous spinal cord stimulation stimulates thedesired motor pool by placing the electrode on top of the vertebrae thatcovers the spinal cord segment of interest. For example, in order tostimulate L2 spinal segment, the electrode is often placed on-top-ofT11-T12 vertebrae segments. However, the applied stimulation current iseasily diverged from the target due to the irregular shape of thevertebrae. This means that the majority of the transcutaneously injectedcurrent flow to other undesired spinal cord segment and thus a highintensity stimulation current is required in order to elicit theresponse of the spinal cord. It is thus unnecessary and ineffective toplace electrodes on top of the desired spinal segment and itscorresponding vertebrae. Instead, the stimulation current can beeffectively injected into the target spinal cord segment through the useof a transcutaneous multi-electrode array that delivers a combination ofdifferent parameters in each electrode. The weighting of the stimulationparameters can be calculated to focus on the segment or interests or toavoid the undesired segment using the described optimization model.

The systems and methods disclosed herein can also be used to stimulatethe spinal cord ganglion with high selectivity, and is applicable todisease associated with movement disorder (e.g., Parkinson's disease,cerebral palsy, and stroke, and even traumatic brain injury) orstimulate sympathetic/non-sympathetic nervous systems.

FIG. 18 shows a schematic diagram of a transcutaneous multi-electrodearray and stimulator system 200 in accordance with the presentdescription. In a preferred embodiment, the electrode array andstimulator system 200 is configured to cover the entire or a portion ofthe spinal cord, including the ganglions and dorsal roots.

The stimulator system 200 comprises an electrode array 202 andcontroller/stimulator 208. The electrode array 202 comprises a flexiblesubstrate 204 housing the array of electrodes 204, and is configured tobe attached onto the patient's back through the 1) use of adhesivehydrogel and/or 2) one or multiple adjustable belts 206 that encirclesthe human chest or abdomen to fix the position of the electrode array.In the array, each electrode is replaceable and can bemounted/dismounted if necessary. One or multiple connectors 206 serve asthe interface to link the electrode array 202 with the transcutaneousstimulator. In a preferred embodiment, the stimulator device 208comprises multi-channel stimulation driver circuits and digitalcontroller (e.g., stimulation shown in system 100 of FIG. 14 )comprising a microprocessor, FPGA, CPU, VPU, and DSP, as well as memory(i.e. volatile or non-volatile memory) to store application programmingin the form of instructions for executing pre-set or real-time receivedstimulation parameters on the microprocessor. The stimulator 208 mayalso comprise a control panel 210 allowing the user to selectivelyconfigure the parameters of each stimulation channel through thestimulator directly. System 250 may also include one or multiplecorresponding connectors (e.g., connectors 206 in FIG. 18 ) from thestimulator link 116 (FIG. 14 ) the stimulator 208 and the electrodearray 202 (FIG. 18 ) for the delivery of the stimulation current. Thestimulator also comprises a power source (not shown) such as arechargeable battery and a wireless module (e.g., Bluetooth and WiFimodules, not shown) for remote control and operation through a mobiledevice (not shown) hold by the patient or clinician.

The electrode array 202 comprises of N×M electrodes 204, where thenumber of the required electrodes varies with different patients and thetargeted anatomical segments. Each electrode has a diameter ranging from0.5 cm to 5 cm and spacing between 0.5 cm to 5 cm. Each electrode 204may be individually replaceable, wherein the deteriorated electrode isremoved and a new electrode is then inserted into the electrode array202.

In one embodiment comprising transcutaneous spinal cord stimulation, thesystem 200 employs a multi-electrode array 202 with M (rows)*N (columns)electrodes. Here M varies from 1 to 100, and N varies from 1 to 100. Theelectrode array 202 generally has a total number varying from 2 to 1000,with diameter ranging between 50 μm and 8 cm. The electrode can have anyshape, such as rectangular, circle, or ring. Note that the electrode 204locations can also be customized, e.g., the electrode 204 arrangementdoes not need to be a rectangular array, but could be in other shapes orarrangements such as circular or occupying a freeform area. The currentapplied to each electrode 204 can be any value provided/calculated bythe optimization algorithm, but should generally be lower than 250 mA,in order to meet the safety requirements. In addition, modulationwaveforms can be provided to support neuromodulation by mixing twowaveforms—e.g., low frequency signals (DC-300 Hz, sinusoidal, squaremono-phasic, square biphasic, triangular) and high frequency signals(high frequency (1 kHz˜40 kHz, sinusoidal, square mono-phasic, squarebiphasic, triangular). In the system of the present description, theelectrode array is able to conduct concurrent recording and stimulationwith the same electrode. The required number, location, and shape of theelectrodes will be determined by the reconstructed spinal cord model andpatient's body shape, multi-physic simulation, and the proposedoptimization method.

In another embodiment comprising transcranial direct/alternating currentstimulation, the disclosed system employs a multi-electrode array 202that total 4˜1000 electrodes 204, with diameter ranging from 0.1 cm˜7cm). For the electrode 204 location, any standard system locations(e.g., international 10-20, 10-10 or 10-5 system), or use customizedlocations can be used. The current applied to each electrode 204 may beany value provided/calculated by the optimization algorithm, but shouldgenerally be lower than 4 mA in order to meet safety requirements duringconstant direct current stimulation. However, one might also use alarger current (e.g., several to tens of mA) with a much shorter pulsewidth (e.g., 50 μs to 10 s of ms). As explained above, the electrodearray is able to conduct concurrent recording and stimulation with thesame electrode 204.

FIG. 19 shows a schematic diagram of a therapeutic spinal cordstimulation system 220 positioned on a patient, integrating thetranscutaneous therapeutic spinal cord stimulator (e.g., stimulator 202of FIG. 18 positioned on back of patient, not shown), vital sign sensors224, and EMG sensors 226, as well as a wireless high-density brainrecording system/stimulation 222 (i.e., EEG, ECoG or LFP/spikerecordings). The therapeutic system 220 is preferably configured as notonly a feedback, but also a feed-forward system. The vital sign sensor224 carried by the patient may be configured to monitors the subject'sheart rate, heart rate variation, EKG, EEG blood pressure/skinimpedance, respiration rate, PPG, SpO2, blood pressure, EGG/ECoG, sweat,etc. are measured to ensure the safety of the stimulation and monitorthe subject's condition. The vital sign sensor 224 is important, asspinal cord stimulation also modulates sympathetic and parasympatheticnerves that travel along with the spinal cord. Once an abnormalphysiological signal is sensed, the stimulation parameters may beadapted or ceased. In another application, the spinal cord stimulationsystem can be used to regulate high blood pressure (BP) and the patent'sBP can be sensed accordingly as a feedback control. In anotherapplication of using spinal cord stimulation to facilitategastrointestinal (GI) motility, placing multiple sensors on theabdominal wall, can record electrogastrogram (EGG) to monitor GImotility.

In one embodiment, The signal from each wireless sensor 224/226 issynchronized and relayed to the stimulator system 200 or externalstorage (not shown) or the electronic health record of thehospital/clinics for feedback or feed-ward control, data analysis andremote patient monitoring.

The EMG sensors and accelerometers 226 are preferably deployed to themuscles of four limbs to sense both flexor and extensor muscle movementsand the patients' postures during the system operation. Thecharacteristics of EMG signals (i.e., the EMG intensity and firingpatterns) and limb postures are used as important biomarkers forevaluating the effectiveness of the stimulation. Not only the earlyresponse of the EMG signal that appears right after stimulation, butalso the middle and late EMG responses appearing several ms or up-to100s of ms after stimulation are used to evaluate the efficacy ofstimulation. It is also important to point out that different motormovements will generally require different stimulation parameters. Forexample, enabling walking and stepping of a paralyzed subject would usedistinct stimulation parameters, and the parameters should be real-timeadapted. The stimulator 200 would thus change its stimulation parametersadaptively based on the sensed EMG signals and patient's posture.

Brain signal recording serves as an important feed-forward controlmechanism in the therapeutic system. Unlike conventional brain-computerinterface or brain-machine-interface technology, recorded brain signalsare directly used to trigger the actuator or served as trigger signal tostimulate the muscles. In the therapeutic system of the presentdescription, the recorded brain signal is used to tune the parameters ofspinal cord stimulation, (i.e. the activated electrode, stimulationfrequencies, intensities, polarities, and pulse width). For instance,when the patient intends to climb the stairs, the recorded brain signalcan be used to tune the stimulation patterns.

In some embodiments, the use of the system of the present descriptionmay be classified into two phases. In the first phase (the trainingphase), the patient/subject uses the device/system/software a couple oftimes per day to discover most efficacious stimulation parameters andstimulation spots, as well as learn of the operation of the system. Inthe second phase, the patient/subject wears the devices/system duringtheir daily activities and sensors will collect the signal of interestsas described previously for further analysis and processing.

FIG. 20 is a schematic block diagram of a wireless sensing andstimulation system/network 250 with a controller 208 configured to beheld by the patient or the clinician. The controller 208 and stimulator200 (or stimulator 100 of FIG. 14 ) are wirelessly coupled to multiplewireless sensors 226, and comprises a mobile device/App with GUI 210designed for remote operation of the stimulator and wireless sensors226. In addition to having the user control the stimulation through thestimulator 200 physically, the wireless controller 208 allowspatient/clinician to configure the stimulator 200 and monitor therecorded physiological signals wirelessly. The recorded physiologicalsignal is sent to either the stimulator 200 or the wireless controller208 (e.g., tablet or like device) for signal analysis and the estimationof optimized stimulation parameters. Once there is a need to update orstop the stimulation parameters, new command is sent either through thestimulator 200 or the mobile device/App. of controller 208. Thecollected data will be eventually uploaded to cloud storage or can besent to the electronic health record used in the hospital.

FIG. 21 is a schematic block diagram of a sensor network in accordancewith the present description. In a preferred embodiment of the presentdescription, the system may include multiple low-power wirelesstransceivers (226a through 226 e) that support duplex data transmission(e.g., Bluetooth) for the sensor network implantation. Instead ofconnecting each wireless sensor to a common central hub, each sensorconnects with each other. For example, 1^(st) sensor 226 a transmits itsrecorded signal to the 2^(nd) sensor 226 b. The 2^(nd) sensor 226 btransmits its own-recorded signal along with signal derived from the1^(st) sensor to 3^(rd) sensor 226 c, and so on. The last sensor #n (226e, which may be one of the EMG sensors 226 or the vital sign sensor 224shown in FIG. 19 ) collects all recorded physiological signals, e.g.,#n−1 sensor 226 d and all previous sensors, and transmits it to thestimulator 100/200 or remote control device 208. As an alternativemethod, each sensor can also connect to a central hop individually andthe transmitted signal contains a time stamp information. The centralhop then collects input from all sensors and synchronize all signalsbased on the timing information provided by each sensor.

Based on existing EM simulation results showing that the skintemperature can be elevated to induce skin burn due to the heataccumulated at the outermost skin layer that has highest resistivity,the duration of each stimulation/therapeutic session should be wellcontrolled to avoid tissue damage and undesired side effects. If apatch-like planar or gel electrode is used for stimulation, eachstimulation session should be limited to tens of minutes (e.g., 1-10mins, 11-20 mins, 21-30 mins, 31-40 mins) to avoid skin burning andpatient discomfort, while the separation between eachstimulation/therapeutic session should be at least >1 minute. If thepenetrating electrode is used, instead, to bypass the high-resistivityskin layer, (see PCT International Application No. PCT/US2016/063886filed on Nov. 28, 2016 and published as WO 2017/091828 A1 on Jun. 1,2017, herein incorporated by reference in its entirety), similarconstraints should be applied to avoid unwanted side effects. Theimpedance of the electrode 204 should also be constantly monitored toensure the robustness and reliability of both electrodes 204 and thestimulator 100/200. The impedance may be measured by measuring the peakvoltage of the electrode overpotential induced by the known stimulationcurrent. Moreover, the impedance between stimulation electrodes may alsobe carefully examined. As the device/system is designed for the dailyuse of patients/subjects, sweating would be an inevitable situation.Impedance between electrodes may therefore be measured to ensure thereis no short circuit current flowing in between.

FIG. 22 shows a schematic diagram of a modulator 230 that modulates thestimulator output 232 and a modulation signal 236 specified by the user.In one embodiment, the modulator 230 comprises a signal mixer or aswitch connected to the output of the stimulator 100/200. The stimulatoroutput waveforms may be selected from low frequency signals (DC-300 Hz,sinusoidal, square mono-phasic, square bi-phasic, triangular signals)and high frequency signals (high frequency (1 kHz-40 kHz, sinusoidal,square mono-phasic, square biphasic, triangular signals), or therecorded electrophysiological signal waveforms (e.g., EMG, EKG, localfield potential waveforms and action potential at its naturalfrequencies)

The user specifies the pulse width and separation of the modulationsignal 234 to determine the number or the length of stimuli allowed topass the mixer/switch per second. The number of allowable stimulationsignals can range from 1 to 20,000 per second based on the efficacy ofthe treatment. The pulse separation can vary from 10 s to 0.001 s whilethe separation and pulse width between each adjacent pulse is adjustableand can be different to create a versatile and flexible stimulationwaveform. This advantageously creates stimuli that comprises of a widerange of frequencies, which might help improve the treatment efficacy,as the physiological signal itself usually includes a wide range offrequencies.

Embodiments of the present technology may be described herein withreference to flowchart illustrations of methods and systems according toembodiments of the technology, and/or procedures, algorithms, steps,operations, formulae, or other computational depictions, which may alsobe implemented as computer program products. In this regard, each blockor step of a flowchart, and combinations of blocks (and/or steps) in aflowchart, as well as any procedure, algorithm, step, operation,formula, or computational depiction can be implemented by various means,such as hardware, firmware, and/or software including one or morecomputer program instructions embodied in computer-readable programcode. As will be appreciated, any such computer program instructions maybe executed by one or more computer processors, including withoutlimitation a general purpose computer or special purpose computer, orother programmable processing apparatus to produce a machine, such thatthe computer program instructions which execute on the computerprocessor(s) or other programmable processing apparatus create means forimplementing the function(s) specified.

Accordingly, blocks of the flowcharts, and procedures, algorithms,steps, operations, formulae, or computational depictions describedherein support combinations of means for performing the specifiedfunction(s), combinations of steps for performing the specifiedfunction(s), and computer program instructions, such as embodied incomputer-readable program code logic means, for performing the specifiedfunction(s). It will also be understood that each block of the flowchartillustrations, as well as any procedures, algorithms, steps, operations,formulae, or computational depictions and combinations thereof describedherein, can be implemented by special purpose hardware-based computersystems which perform the specified function(s) or step(s), orcombinations of special purpose hardware and computer-readable programcode.

Furthermore, these computer program instructions, such as embodied incomputer-readable program code, may also be stored in one or morecomputer-readable memory or memory devices that can direct a computerprocessor or other programmable processing apparatus to function in aparticular manner, such that the instructions stored in thecomputer-readable memory or memory devices produce an article ofmanufacture including instruction means which implement the functionspecified in the block(s) of the flowchart(s). The computer programinstructions may also be executed by a computer processor or otherprogrammable processing apparatus to cause a series of operational stepsto be performed on the computer processor or other programmableprocessing apparatus to produce a computer-implemented process such thatthe instructions which execute on the computer processor or otherprogrammable processing apparatus provide steps for implementing thefunctions specified in the block(s) of the flowchart(s), procedure (s)algorithm(s), step(s), operation(s), formula(e), or computationaldepiction(s).

It will further be appreciated that the terms “programming” or “programexecutable” as used herein refer to one or more instructions that can beexecuted by one or more computer processors to perform one or morefunctions as described herein. The instructions can be embodied insoftware, in firmware, or in a combination of software and firmware. Theinstructions can be stored local to the device in non-transitory media,or can be stored remotely such as on a server, or all or a portion ofthe instructions can be stored locally and remotely. Instructions storedremotely can be downloaded (pushed) to the device by user initiation, orautomatically based on one or more factors.

It will further be appreciated that as used herein, that the termsprocessor, hardware processor, computer processor, central processingunit (CPU), and computer are used synonymously to denote a devicecapable of executing the instructions and communicating withinput/output interfaces and/or peripheral devices, and that the termsprocessor, hardware processor, computer processor, CPU, and computer areintended to encompass single or multiple devices, single core andmulticore devices, and variations thereof.

From the description herein, it will be appreciated that the presentdisclosure encompasses multiple embodiments which include, but are notlimited to, the following:

1. A method for optimized transcutaneous stimulation of a targettreatment region, the method comprising: generating a 3D model of thetarget treatment region as a function of an MRI or CT image of thetarget treatment region, a stimulation electrode array to be applied tothe target treatment region, and calculation of a lead field matrixassociated with the stimulation electrode array and the target treatmentregion; determining a safety limit for transcutaneous stimulation of thetarget treatment region; and generating an optimization model forstimulation parameters that provide both high intensity and focalaccuracy within safety the safety limit by applying a parameter thatassigns a weight to a directional intensity and focality associated withsaid transcutaneous stimulation; wherein said method is performed by aprocessor executing instructions stored on a non-transitory medium.

2. The method, apparatus or system of any preceding or subsequentembodiment, wherein generating a 3D model of the target treatment regioncomprises: acquiring the MRI or CT image of the target treatment region;segmenting the MRI or CT image into different tissues according to greylevels within the MRI or CT image; generating a target tissue model ofthe target treatment region for each of the different tissues;constructing an electrode model based on the stimulation electrode arrayand the constructed target tissue model; discretizing the target tissuemodel and electrode model into a large number of voxels to form a finiteelement model; and calculating a lead field matrix of the finite elementmodel.

3. The method, apparatus or system of any preceding or subsequentembodiment, wherein the target treatment comprises the spine foroptimization of transcutaneous spinal cord stimulation (tSCS).

4. The method, apparatus or system of any preceding or subsequentembodiment, wherein the target treatment comprises the brain fortranscranial current stimulation (tCS).

5. The method, apparatus or system of any preceding or subsequentembodiment, wherein the target treatment comprises a region of the torsofor internal organ stimulation.

6. The method, apparatus or system of any preceding or subsequentembodiment, wherein determining a safety limit is performed according tothe function:

$\quad\left\{ {\begin{matrix}{{{x_{i}} \leq I_{\max}},{i = 1},\ldots\;,N} \\{{\sum\limits_{i = 1}^{N}{x_{i}}} \leq {2*I_{total}}} \\{{\sum\limits_{i = 1}^{N}x_{i}} = 0} \\{{Intensity}_{avoid} \leq {\frac{1}{ratio}{Intensity}_{target}}}\end{matrix};} \right.$

wherein I_(max) represents a maximum current at each electrode,I_(total) denotes a maximum total injected current, and ratio representsan intensity ratio between the target treatment region and an avoidanceregion.

7. The method, apparatus or system of any preceding or subsequentembodiment, wherein generating an optimization model is performedaccording to:

$\begin{matrix}{{x_{proposed} = {{\arg\mspace{11mu}{\min\limits_{x}{\frac{1}{w}{{Kx}}_{2}^{2}}}} - {\lambda*e_{0}^{T}{Cx}}}}{{subject}\mspace{14mu}{to}\mspace{14mu}\left\{ \begin{matrix}{{{x_{i}} \leq I_{\max}},{i = 1},\ldots\;,N} \\{{\sum\limits_{i = 1}^{N}{x_{i}}} \leq {2*I_{total}}} \\{{\sum\limits_{i = 1}^{N}x_{i}} = 0} \\{{Intensity}_{avoid} \leq {\frac{1}{ratio}{Intensity}_{target}}}\end{matrix} \right.}} & \;\end{matrix}$wherein w is a constant equal to a ratio between a total number ofvoxels and a number of targeted voxels; wherein

$\arg\mspace{11mu}{\min\limits_{x}{\frac{1}{w}{{Kx}}_{2}^{2}}}$applies the focality and λ*e₀ ^(T) Cx applies to the directionalintensity; and wherein λ is the parameter that assigns a weight to adirectional intensity and focality.

8. The method, apparatus or system of any preceding or subsequentembodiment: wherein by varying λ, an upper bound of focality anddirectional intensity can be estimated; and wherein directionalintensity increases as λ increases and focality increases as λdecreases.

9. The method, apparatus or system of any preceding or subsequentembodiment, further comprising: performing EEG/EMG inverse image guidedoptimal stimulation to the target treatment region with the stimulationelectrode array.

10. The method, apparatus or system of any preceding or subsequentembodiment, wherein the stimulation parameters provide direct currentstimulation for focalized stimulation of the target tissue region at anyorientation with high precision.

11. The method, apparatus or system of any preceding or subsequentembodiment, wherein multiple targets within the target region can besimultaneously stimulated or designated regions within the target tissueregion can be avoided.

12. An apparatus for optimized transcutaneous stimulation of a targettreatment region, the apparatus comprising: (a) a processor; and (b) anon-transitory memory storing instructions executable by the processor;(c) wherein said instructions, when executed by the processor, performsteps comprising: (i) generating a 3D model of the target treatmentregion as a function of an MRI or CT image of the target treatmentregion, a stimulation electrode array to be applied to the targettreatment region, and calculation of a lead field matrix associated withthe stimulation electrode array and the target treatment region; (ii)generating an optimization model for stimulation parameters that provideboth high intensity and focal accuracy within safety the safety limit byapplying a parameter that assigns a weight to a directional intensityand focality associated with said transcutaneous stimulation; and (iii)determining a safety limit for transcutaneous stimulation of the targettreatment region.

13. The method, apparatus or system of any preceding or subsequentembodiment, wherein generating a 3D model of the target treatment regioncomprises: acquiring the MRI or CT image of the target treatment region;segmenting the MRI or CT image into different tissues according to greylevels within the MRI or CT image; generating a target tissue model ofthe target treatment region for each of the different tissues;constructing an electrode model based on the stimulation electrodearray; discretizing the target tissue model and electrode model into alarge number of voxels to form a finite element model; and calculating alead field matrix of the finite element model.

14. The method, apparatus or system of any preceding or subsequentembodiment, wherein the target treatment comprises the spine foroptimization of transcutaneous spinal cord stimulation (tSCS).

15. The method, apparatus or system of any preceding or subsequentembodiment, wherein the target treatment comprises the brain fortranscranial current stimulation (tCS).

16. The method, apparatus or system of any preceding or subsequentembodiment, wherein the target treatment comprises a region of the torsofor internal organ stimulation.

17. The method, apparatus or system of any preceding or subsequentembodiment, wherein determining a safety limit is performed according tothe function

$\left\{ {\begin{matrix}{{{x_{i}} \leq I_{\max}},{i = 1},\ldots\;,N} \\{{\sum\limits_{i = 1}^{N}{x_{i}}} \leq {2*I_{total}}} \\{{\sum\limits_{i = 1}^{N}x_{i}} = 0} \\{{Intensity}_{avoid} \leq {\frac{1}{ratio}{Intensity}_{target}}}\end{matrix};} \right.$wherein I_(max) represents a maximum current at each electrode,I_(total) denotes a maximum total injected current, and ratio representsan intensity ratio between the target treatment region and an avoidanceregion.

18. The method, apparatus or system of any preceding or subsequentembodiment, wherein generating an optimization model is performedaccording to:

$\begin{matrix}{{x_{proposed} = {{\arg\mspace{11mu}{\min\limits_{x}{\frac{1}{w}{{Kx}}_{2}^{2}}}} - {\lambda*e_{0}^{T}{Cx}}}}{{subject}\mspace{14mu}{to}\mspace{14mu}\left\{ \begin{matrix}{{{x_{i}} \leq I_{\max}},{i = 1},\ldots\;,N} \\{{\sum\limits_{i = 1}^{N}{x_{i}}} \leq {2*I_{total}}} \\{{\sum\limits_{i = 1}^{N}x_{i}} = 0} \\{{Intensity}_{avoid} \leq {\frac{1}{ratio}{Intensity}_{target}}}\end{matrix} \right.}} & \;\end{matrix}$wherein w is a constant equal to a ratio between a total number ofvoxels and a number of targeted voxels; wherein

$\arg\mspace{11mu}{\min\limits_{x}{\frac{1}{w}{{Kx}}_{2}^{2}}}$applies the focality and λ*e₀ ^(T)Cx applies to the directionalintensity; and wherein λ is the parameter that assigns a weight to adirectional intensity and focality.

19. The method, apparatus or system of any preceding or subsequentembodiment: wherein by varying λ, an upper bound of focality anddirectional intensity can be estimated; and wherein directionalintensity increases as λ increases and focality increases as λdecreases.

20. The method, apparatus or system of any preceding or subsequentembodiment, wherein said instructions when executed by the processorfurther perform steps comprising: performing EEG/EMG inverse imageguided optimal stimulation to the target treatment region with thestimulation electrode array.

21. The method, apparatus or system of any preceding or subsequentembodiment, wherein the stimulation parameters provide direct currentstimulation or alternative current stimulation at a short period of timefor focalized stimulation of the target tissue region at any orientationwith high precision.

22. The method, apparatus or system of any preceding or subsequentembodiment, wherein multiple targets within the target region can besimultaneously stimulated or a designated region within the targettissue region can be avoided.

23. The method, apparatus or system of any preceding or subsequentembodiment, further comprising: a multi-electrode array comprising anetwork of high-density electrodes; and a stimulator comprising highvoltage CMOS/BJT devices that are employed as current sources; thestimulator comprising a global controller and a plurality of localcontrollers, both the global controller and plurality of localcontrollers being configured to accept data packets for a stimulationconfiguration over a single input.

24. The method, apparatus or system of any preceding or subsequentembodiment, the stimulator further configured to simultaneously anddigitally program each electrode channel simultaneously.

25. The method, apparatus or system of any preceding or subsequentembodiment, the stimulator further configured to evoke a 2-D electricfield pattern in-vitro by selectively choosing the stimulation channels,intensities, polarities, and return channel from selected electrodeswithin the electrode array.

26. A system for optimized transcutaneous stimulation of a targettreatment region, the system comprising: (a) a multi-electrode arraycomprising a network of high-density electrodes; (b) a stimulatorconfigured to control operation of the multi-electrode array; (c) aprocessor; and (d) a non-transitory memory storing instructionsexecutable by the processor; (e) wherein said instructions, whenexecuted by the processor, perform steps comprising; (i) generating a 3Dmodel of the target treatment region as a function of an MRI or CT imageof the target treatment region, a stimulation electrode array to beapplied to the target treatment region, and calculation of a lead fieldmatrix associated with the stimulation electrode array and the targettreatment region; (ii) generating an optimization model for stimulationparameters that provide both high intensity and focal accuracy withinsafety the safety limit by applying a parameter that assigns a weight toa directional intensity and focality associated with said transcutaneousstimulation; and (iii) determining a safety limit for transcutaneousstimulation of the target treatment region; and (iv) performing EEG/EMGinverse image guided optimal stimulation to the target treatment regionwith the stimulation electrode array.

27. The method, apparatus or system of any preceding or subsequentembodiment, wherein generating a 3D model of the target treatment regioncomprises: acquiring the MRI or CT image of the target treatment region;segmenting the MRI or CT image into different tissues according to greylevels within the MRI or CT image; generating a target tissue model ofthe target treatment region for each of the different tissues;constructing an electrode model based on the stimulation electrodearray; discretizing the target tissue model and electrode model into alarge number of voxels to form a finite element model; and calculating alead field matrix of the finite element model.

28. The method, apparatus or system of any preceding or subsequentembodiment, wherein the target treatment comprises the spine foroptimization of transcutaneous spinal cord stimulation (tSCS).

29. The method, apparatus or system of any preceding or subsequentembodiment, wherein the target treatment comprises the brain fortranscranial current stimulation (tCS).

30. The method, apparatus or system of any preceding or subsequentembodiment, wherein the target treatment comprises the a region of thetorso for internal organ stimulation.

31. The method, apparatus or system of any preceding or subsequentembodiment, wherein determining a safety limit is performed according tothe function

$\quad\left\{ {\begin{matrix}{{{x_{i}} \leq I_{\max}},{i = 1},\ldots\;,N} \\{{\sum\limits_{i = 1}^{N}{x_{i}}} \leq {2*I_{total}}} \\{{\sum\limits_{i = 1}^{N}x_{i}} = 0} \\{{Intensity}_{avoid} \leq {\frac{1}{ratio}{Intensity}_{target}}}\end{matrix};} \right.$wherein I_(max) represents a maximum current at each electrode,I_(total) denotes a maximum total injected current, and ratio representsan intensity ratio between the target treatment region and an avoidanceregion.

32. The method, apparatus or system of any preceding or subsequentembodiment, wherein generating an optimization model is performedaccording to:

$\begin{matrix}{{x_{proposed} = {{\arg\mspace{11mu}{\min\limits_{x}{\frac{1}{w}{{Kx}}_{2}^{2}}}} - {\lambda*e_{0}^{T}{Cx}}}}{{subject}\mspace{14mu}{to}\mspace{14mu}\left\{ \begin{matrix}{{{x_{i}} \leq I_{\max}},{i = 1},\ldots\;,N} \\{{\sum\limits_{i = 1}^{N}{x_{i}}} \leq {2*I_{total}}} \\{{\sum\limits_{i = 1}^{N}x_{i}} = 0} \\{{Intensity}_{avoid} \leq {\frac{1}{ratio}{Intensity}_{target}}}\end{matrix} \right.}} & \;\end{matrix}$wherein w is a constant equal to a ratio between a total number ofvoxels and a number of targeted voxels; wherein

$\arg\mspace{11mu}{\min\limits_{x}{\frac{1}{w}{{Kx}}_{2}^{2}}}$applies the focality and λ*e₀ ^(T)Cx applies to the directionalintensity; and wherein λ is the parameter that assigns a weight to adirectional intensity and focality.

33. The method, apparatus or system of any preceding or subsequentembodiment: wherein by varying λ, an upper bound of focality anddirectional intensity can be estimated; and wherein directionalintensity increases as λ increases and focality increases as λdecreases.

34. The method, apparatus or system of any preceding or subsequentembodiment, wherein the stimulation parameters provide direct currentstimulation for focalized stimulation of the target tissue region at anyorientation with high precision.

35. The method, apparatus or system of any preceding or subsequentembodiment, wherein multiple targets within the target region can besimultaneously stimulated or a designated region within the targettissue region can be avoided.

36. The method, apparatus or system of any preceding or subsequentembodiment, wherein the stimulator comprises a global controller and aplurality of local controllers, both the global controller and pluralityof local controllers being configured to accept data packets for astimulation configuration over a single input.

37. The method, apparatus or system of any preceding or subsequentembodiment, the stimulator further configured to simultaneously anddigitally program each electrode channel simultaneously.

38. The method, apparatus or system of any preceding or subsequentembodiment, the stimulator further configured to evoke a 2-D or 3-Delectric field pattern in-vitro by selectively choosing the stimulationchannels, intensities, polarities, and return channel from selectedelectrodes within the electrode array.

39. The method, apparatus or system of any preceding or subsequentembodiment, wherein the return electrode is selected from a single or asubset of the electrodes in the array or by placing a return on anopposing side of the target tissue region.

40. The apparatus of any preceding or subsequent embodiment, thestimulator configured for generating a stimulation waveform by mixing atleast waveforms from among low frequency signals, high frequency signalsor the electrophysiological signal waveforms.

41. The method, apparatus or system of any preceding or subsequentembodiment, wherein the size of the electrodes ranges from 0.5 cm to 5cm and the spacing between each electrode ranges from 0.5 cm to 5 cm.

42. The apparatus or method comprising using optimization techniques toimprove the focal accuracy of transcutaneous spinal cord stimulation(tSCS), transcranial current stimulation (tCS), or internal organstimulation.

43. The method, apparatus or system of any preceding or subsequentembodiment, wherein high intensity and focality of a target is achieved.

44. A method for optimizing stimulation parameters for transcutaneousspinal cord stimulation (tSCS), the method comprising: (a) determining asafety limit for transcutaneous spinal cord stimulation according to

$\quad\left\{ \begin{matrix}{{{x_{i}} \leq I_{\max}},{i = 1},\ldots\;,N} & (1) \\{{\sum\limits_{i = 1}^{N}{x_{i}}} \leq {2*I_{total}}} & (2) \\{{\sum\limits_{i = 1}^{N}x_{i}} = 0} & (3) \\{{Intensity}_{avoid} \leq {\frac{1}{ratio}{Intensity}_{target}}} & (4)\end{matrix} \right.$where I_(max) represents the maximum current at each electrode,I_(total) denotes the maximum total current injected to the body, andratio represents the intensity ratio between the target and avoidanceregion; (b) developing an optimization model that always provides afeasible solution, which provides both high intensity and focal accuracywithin safety constraints according to

$\begin{matrix}{{x_{proposed} = {{\arg\mspace{11mu}{\min\limits_{x}{\frac{1}{w}{{Kx}}_{2}^{2}}}} - {\lambda*e_{0}^{T}{Cx}}}}{{subject}\mspace{14mu}{to}\mspace{14mu}\left\{ \begin{matrix}{{{x_{i}} \leq I_{\max}},{i = 1},\ldots\;,N} \\{{\sum\limits_{i = 1}^{N}{x_{i}}} \leq {2*I_{total}}} \\{{\sum\limits_{i = 1}^{N}x_{i}} = 0} \\{{Intensity}_{avoid} \leq {\frac{1}{ratio}{Intensity}_{target}}}\end{matrix} \right.}} & \;\end{matrix}$where the constant w is equal to the ratio between the total number ofvoxels and the number of targeted voxels, where the first term is thefocality term and the second term is the intensity on the desireddirection, and where the parameter λ balances these two objectives andcontrols the relative importance of the focality and directionalintensity; (c) performing EEG/EMG Inverse Image Guided OptimalStimulation applied to transcutaneous spinal cord stimulation; and (d)using the foregoing to design optimal parameters for transcutaneousspinal cord stimulation to achieve focalized stimulation wherein anytarget location can be stimulated with any orientations with highprecision, and wherein single or multiple targets can be stimulated andwherein designated regions can be avoided.

45. The method, apparatus or system of any preceding or subsequentembodiment, wherein the parameter λ controls the relative importancebetween the intensity and the focality.

46. The method, apparatus or system of any preceding or subsequentembodiment, wherein by changing λ an upper bound of focality andintensity is estimated.

47. The method, apparatus or system of any preceding or subsequentembodiment: wherein when λ is very large, the best intensity results;and wherein when λ is very small, the best focality results.

48. A method for optimizing stimulation parameters for transcranialcurrent stimulation (tCS), the method comprising: (a) determining asafety limit for transcranial direct current stimulation according to

$\quad\left\{ \begin{matrix}{{{x_{i}} \leq I_{\max}},{i = 1},\ldots\;,N} & (1) \\{{\sum\limits_{i = 1}^{N}{x_{i}}} \leq {2*I_{total}}} & (2) \\{{\sum\limits_{i = 1}^{N}x_{i}} = 0} & (3) \\{{Intensity}_{avoid} \leq {\frac{1}{ratio}{Intensity}_{target}}} & (4)\end{matrix} \right.$where I_(max) represents the maximum current at each electrode,I_(total) denotes the maximum total current injected to the body, andratio represents the intensity ratio between the target and avoidanceregion; (b) developing an optimization model that always provides afeasible solution, which provides both high intensity and focal accuracywithin safety constraints according to

$x_{proposed} = {{\arg{\;\;}{\min\limits_{x}\;{\frac{1}{w}{{Kx}}_{2}^{2}}}} - {\lambda*e_{0}^{T}{Cx}}}$${subject}\mspace{14mu}{to}\mspace{14mu}\left\{ \begin{matrix}{{{x_{i}} \leq I_{\max}},{i = 1},\ldots\;,N} \\{{\sum\limits_{i = 1}^{N}{x_{i}}} \leq {2*I_{total}}} \\{{\sum\limits_{i = 1}^{N}x_{i}} = 0} \\{{Intensity}_{avoid} \leq {\frac{1}{ratio}{Intensity}_{target}}}\end{matrix} \right.$where the constant w is equal to the ratio between the total number ofvoxels and the number of targeted voxels, where the first term is thefocality term and the second term is the intensity on the desireddirection, and where the parameter λ balances these two objectives andcontrols the relative importance of the focality and directionalintensity; (c) performing EEG/EMG Inverse Image Guided OptimalStimulation applied to transcranial current stimulation; and (d) usingthe foregoing to design optimal parameters for transcranial directcurrent stimulation to achieve focalized stimulation wherein any targetlocation can be stimulated with any orientations with high precision,and wherein single or multiple targets can be stimulated and whereindesignated regions can be avoided.

49. The method, apparatus or system of any preceding or subsequentembodiment, wherein the parameter λ controls the relative importancebetween the intensity and the focality.

50. The method, apparatus or system of any preceding or subsequentembodiment, wherein by changing λ an upper bound of focality andintensity is estimated.

51. The method, apparatus or system of any preceding or subsequentembodiment: wherein when λ is very large, the best intensity results;and wherein when λ is very small, the best focality results.

52. An apparatus for optimizing stimulation parameters fortranscutaneous spinal cord stimulation (tSCS), the apparatus comprising:(a) a multi-electrode array; (b) a display device; (c) a hardwareprocessor; (d) a non-transitory memory storing instructions executableby the hardware processor; (e) said instructions, when executed,performing steps comprising: (i) determining a safety limit fortranscutaneous spinal cord stimulation according to

$\quad\left\{ \begin{matrix}{{{x_{i}} \leq I_{\max}},{i = 1},\ldots\;,N} & (1) \\{{\sum\limits_{i = 1}^{N}{x_{i}}} \leq {2*I_{total}}} & (2) \\{{\sum\limits_{i = 1}^{N}x_{i}} = 0} & (3) \\{{Intensity}_{avoid} \leq {\frac{1}{ratio}{Intensity}_{target}}} & (4)\end{matrix} \right.$where I_(max) represents the maximum current at each electrode,I_(total) denotes the maximum total current injected to the body, andratio represents the intensity ratio between the target and avoidanceregion; (ii) developing an optimization model that always provides afeasible solution, which provides both high intensity and focal accuracywithin safety constraints according to

$x_{proposed} = {{\arg{\;\;}{\min\limits_{x}\;{\frac{1}{w}{{Kx}}_{2}^{2}}}} - {\lambda*e_{0}^{T}{Cx}}}$${subject}\mspace{14mu}{to}\mspace{14mu}\left\{ \begin{matrix}{{{x_{i}} \leq I_{\max}},{i = 1},\ldots\;,N} \\{{\sum\limits_{i = 1}^{N}{x_{i}}} \leq {2*I_{total}}} \\{{\sum\limits_{i = 1}^{N}x_{i}} = 0} \\{{Intensity}_{avoid} \leq {\frac{1}{ratio}{Intensity}_{target}}}\end{matrix} \right.$where the constant w is equal to the ratio between the total number ofvoxels and the number of targeted voxels, where the first term is thefocality term and the second term is the intensity on the desireddirection, and where the parameter λ balances these two objectives andcontrols the relative importance of the focality and directionalintensity; (iii) performing EEG/EMG Inverse Image Guided OptimalStimulation applied to transcutaneous spinal cord stimulation; (iv)using the foregoing to design optimal parameters for transcutaneousspinal cord stimulation to achieve focalized stimulation wherein anytarget location can be stimulated with any orientations with highprecision, and wherein single or multiple targets can be stimulated andwherein designated regions can be avoided; and (v) displayingstimulation information on the display device.

53. The method, apparatus or system of any preceding or subsequentembodiment, wherein the parameter λ controls the relative importancebetween the intensity and the focality.

54. The method, apparatus or system of any preceding or subsequentembodiment, wherein by changing λ an upper bound of focality andintensity is estimated.

55. The method, apparatus or system of any preceding or subsequentembodiment: wherein when λ is very large, the best intensity results;and wherein when λ is very small, the best focality results.

56. An apparatus for optimizing stimulation parameters for transcranialcurrent stimulation (tCS), the apparatus comprising: (a) amulti-electrode array; (b) a display device; (c) a hardware processor;(d) a non-transitory memory storing instructions executable by thehardware processor; (e) said instructions, when executed, performingsteps comprising: (i) determining a safety limit for transcranial directcurrent stimulation according to

$\quad\left\{ \begin{matrix}{{{x_{i}} \leq I_{\max}},{i = 1},\ldots\;,N} & (1) \\{{\sum\limits_{i = 1}^{N}{x_{i}}} \leq {2*I_{total}}} & (2) \\{{\sum\limits_{i = 1}^{N}x_{i}} = 0} & (3) \\{{Intensity}_{avoid} \leq {\frac{1}{ratio}{Intensity}_{target}}} & (4)\end{matrix} \right.$where I_(max) represents the maximum current at each electrode,I_(total) denotes the maximum total current injected to the body, andratio represents the intensity ratio between the target and avoidanceregion; (ii) developing an optimization model that always provides afeasible solution, which provides both high intensity and focal accuracywithin safety constraints according to

$x_{proposed} = {{\arg\;{\underset{x}{\;\min}\;{\frac{1}{w}{{Kx}}_{2}^{2}}}} - {\lambda*e_{0}^{T}{Cx}}}$${subject}\mspace{14mu}{to}\mspace{14mu}\left\{ \begin{matrix}{{{x_{i}} \leq I_{\max}},{i = 1},\ldots\;,N} \\{{\sum\limits_{i = 1}^{N}{x_{i}}} \leq {2*I_{total}}} \\{{\sum\limits_{i = 1}^{N}x_{i}} = 0} \\{{Intensity}_{avoid} \leq {\frac{1}{ratio}{Intensity}_{target}}}\end{matrix} \right.$where the constant w is equal to the ratio between the total number ofvoxels and the number of targeted voxels, where the first term is thefocality term and the second term is the intensity on the desireddirection, and where the parameter λ balances these two objectives andcontrols the relative importance of the focality and directionalintensity; (iii) performing EEG/EMG Inverse Image Guided OptimalStimulation applied to transcranial current stimulation; (iv) using theforegoing to design optimal parameters for transcranial direct currentstimulation to achieve focalized stimulation wherein any target locationcan be stimulated with any orientations with high precision, and whereinsingle or multiple targets can be stimulated and wherein designatedregions can be avoided; and (v) displaying stimulation information onthe display device.

57. The method, apparatus or system of any preceding or subsequentembodiment, wherein the parameter λ controls the relative importancebetween the intensity and the focality.

58. The method, apparatus or system of any preceding or subsequentembodiment, wherein by changing λ an upper bound of focality andintensity is estimated.

59. The method, apparatus or system of any preceding or subsequentembodiment: wherein when λ is very large, the best intensity results;and wherein when λ is very small, the best focality results.

60. An apparatus for non-invasive spinal cord stimulation for motilitydisorder, the apparatus comprising: a wireless controller; a stimulator;and a plurality of wireless sensors; wherein said controller, saidstimulator, and said wireless sensors are configured for: allowing auser to remotely operate the stimulator; allowing the user or aclinician to configure the stimulator and monitor the recordedphysiological signals wirelessly; and receiving a recorded physiologicalsignal that is sent to either the stimulator or the wireless controldevice for signal analysis and estimation of optimized stimulationparameters.

61. The method, apparatus or system of any preceding or subsequentembodiment: wherein the user first takes MRI to derive his/her MRIimage; wherein a mechanical model based on MRI images is built,including cerebrospinal fluid (CSF), spinal cord (gray and whitematters), vertebrae, muscle, skin (e.g., stratum germinativum andstratum corneum), as well as fat tissues; and wherein each physicallayer is assigned a corresponding electrical and thermal property (e.g.,conductivity and permittivity) for simulation.

62. The method, apparatus or system of any preceding or subsequentembodiment: wherein the stimulation parameters are adapted based on thephysiological signals captured by sensors and the postures/motions thepatient/subject intends to perform and the stimulation is ceased onceabnormal vital signs are captured by the wireless sensors.

63. The method, apparatus or system of any preceding or subsequentembodiment: wherein the electrode for spinal cord stimulation is aplanar electrode or needle electrode that penetrates dead skin todeliver electrical stimuli.

64. The method, apparatus or system of any preceding or subsequentembodiment, wherein the spinal cord stimulation system is used tomodulate center/peripheral/autonomic nervous system, such as but notlimited to motor function, brain state, GI motility, blood pressure,heart rate, and respiration rate

As used herein, the singular terms “a,” “an,” and “the” may includeplural referents unless the context clearly dictates otherwise.Reference to an object in the singular is not intended to mean “one andonly one” unless explicitly so stated, but rather “one or more.”

As used herein, the term “set” refers to a collection of one or moreobjects. Thus, for example, a set of objects can include a single objector multiple objects.

As used herein, the terms “substantially” and “about” are used todescribe and account for small variations. When used in conjunction withan event or circumstance, the terms can refer to instances in which theevent or circumstance occurs precisely as well as instances in which theevent or circumstance occurs to a close approximation. When used inconjunction with a numerical value, the terms can refer to a range ofvariation of less than or equal to ±10% of that numerical value, such asless than or equal to ±5%, less than or equal to ±4%, less than or equalto ±3%, less than or equal to ±2%, less than or equal to ±1%, less thanor equal to ±0.5%, less than or equal to ±0.1%, or less than or equal to±0.05%. For example, “substantially” aligned can refer to a range ofangular variation of less than or equal to ±10°, such as less than orequal to ±5°, less than or equal to ±4°, less than or equal to ±3°, lessthan or equal to ±2°, less than or equal to ±1°, less than or equal to±0.5°, less than or equal to ±0.1°, or less than or equal to ±0.05°.

Additionally, amounts, ratios, and other numerical values may sometimesbe presented herein in a range format. It is to be understood that suchrange format is used for convenience and brevity and should beunderstood flexibly to include numerical values explicitly specified aslimits of a range, but also to include all individual numerical valuesor sub-ranges encompassed within that range as if each numerical valueand sub-range is explicitly specified. For example, a ratio in the rangeof about 1 to about 200 should be understood to include the explicitlyrecited limits of about 1 and about 200, but also to include individualratios such as about 2, about 3, and about 4, and sub-ranges such asabout 10 to about 50, about 20 to about 100, and so forth.

Although the description herein contains many details, these should notbe construed as limiting the scope of the disclosure but as merelyproviding illustrations of some of the presently preferred embodiments.Therefore, it will be appreciated that the scope of the disclosure fullyencompasses other embodiments which may become obvious to those skilledin the art.

All structural and functional equivalents to the elements of thedisclosed embodiments that are known to those of ordinary skill in theart are expressly incorporated herein by reference and are intended tobe encompassed by the present claims. Furthermore, no element,component, or method step in the present disclosure is intended to bededicated to the public regardless of whether the element, component, ormethod step is explicitly recited in the claims. No claim element hereinis to be construed as a “means plus function” element unless the elementis expressly recited using the phrase “means for”. No claim elementherein is to be construed as a “step plus function” element unless theelement is expressly recited using the phrase “step for”.

TABLE 1 Conductivity Values Used For The 3D Spinal Cord Model ElectricalSC 0.0018 conductivity σ (S) SG 0.15 Dermis 1.62 Fat 0.02383 Abdomen0.0414 Vertebrae 0.0285 Muscle 0.08 (transversal)/0.5 (longitudinal) CSF1.7 Gray matter 0.023 White matter 0.083 (transversal)/0.6(longitudinal)

TABLE 2 Stimulator SoC Summary Max # of stimulation Channels 256 Outcurrent per channel (effect # of channels: 256) −0.5 to 0.5 mA Outcurrent per channel (effect # of channels: 64) −2 to 2 mA Stimulatorcompliance voltage ±15 Stimulator SoC size 5.7 × 6.6 mm² Static powerconsumption 12 mW

What is claimed is:
 1. A method comprising: receiving, by an external device comprising a processor, high-resolution structural images of a target treatment region of a patient, wherein a stimulation electrode array is placed within the target treatment region; generating, by the external device, a 3D model based on high-resolution structural images of the target treatment region with the stimulation electrode array positioned within the target treatment region; calculating, by the external device, based on the 3D model, a lead field matrix associated with the stimulation electrode array and the target treatment region; determining, by the external device, a safety limit for transcutaneous stimulation of the target treatment region by the stimulation electrode array; generating, by the external device, a set of parameters for a stimulation to be delivered by the stimulation electrode array based on the lead field matrix and the safety limit; and communicating by the external device to the stimulation electrode array, the set of stimulation parameters, wherein the stimulation electrode array is configured to deliver a stimulation with high intensity and focal accuracy of electrical stimulation within the safety limit based on the set of stimulation parameters to the target treatment region, wherein the set of stimulation parameters assigns a relative weight to directional intensity and focality of stimulation at the target treatment region.
 2. The method of claim 1, wherein generating the 3D model further comprises: segmenting the high-resolution structural images into different tissues according to grey levels; generating a target tissue model of the target treatment region for each of the different tissues; constructing an electrode model based on the stimulation electrode array and the constructed target tissue model; and discretizing the target tissue model and electrode model into a large number of voxels to form a finite element model, wherein the lead field matrix is calculated based on the finite element model.
 3. The method of claim 1, wherein the target treatment region comprises a spinal cord or spinal nerves.
 4. The method of claim 1, wherein the target treatment region comprises a brain structure.
 5. The method of claim 1, wherein the target treatment region comprises an internal organ.
 6. The method of claim 1, wherein determining the safety limit for transcutaneous stimulation of the target treatment region by the stimulation electrode array further comprises determining a ratio between current to be applied to the target treatment region and to an avoidance region.
 7. The method of claim 1, further comprising: utilizing, by the external device, a stimulation optimization model to provide a feasible solution that includes both the high intensity and focal accuracy of electrical stimulation within the safety limit; and wherein the stimulation optimization model is based on the focality of stimulation defining a ratio between the total number of voxels and the number of targeted voxels and the directional intensity of stimulation in the desired direction, wherein the focality and directional intensity of stimulation are balanced by assigning the relative weight terms.
 8. The method of claim 1, further comprising: measuring, by one or more sensing devices, one or more physiological signals, which includes at least one of EEG, EMG, EKG, EcoG, LFP, spike recording, accelerometer, PPG, SpO₂, respiration, or skin impedance signals; adjusting, by the external device, one or more changes to the set of stimulation parameters in response to the measured physiological signals; and generating, by the external device, an updated set of optimal stimulation parameters.
 9. The method of claim 1, wherein the high-resolution structural images comprise MRI or CT images.
 10. The method of claim 1, wherein the stimulation having the set of stimulation parameters provides direct current stimulation for focalized stimulation of the target treatment region at any orientation with high precision.
 11. The method of claim 1, further comprising designating, simultaneously via the external device, at least one target within the target treatment region as a simulation target and at least one region within the target treatment region as a region to be avoided.
 12. An apparatus comprising: a stimulation electrode array configured to be placed within a target treatment region of a patient; and an external device comprising: a processor; and a non-transitory memory storing instructions executable by the processor to perform steps comprising: receiving high-resolution structural images of the target treatment region; generating a 3D model based on the high-resolution structural images of the target treatment region with the stimulation electrode array placed within the target treatment region; calculating, based on the 3D model, a lead field matrix associated with the stimulation electrode array and the target treatment region; determining a safety limit for transcutaneous stimulation of the target treatment region by the stimulation electrode array; generating a set of stimulation parameters to be delivered by the stimulation electrode array based on the lead field matrix and the safety limit; and communicating to the stimulation electrode array the set of stimulation parameters, wherein the stimulation electrode array is configured to deliver a stimulation with high intensity and focal accuracy of electrical stimulation within the safety limit based on the set of stimulation parameters, where the set of stimulation parameters assigns a relative weight to a directional intensity and a focality of stimulation of the target treatment region.
 13. The apparatus of claim 12, wherein the generating a 3D model further comprises: segmenting the high-resolution structural images into different tissues according to grey levels; generating a target tissue model of the target treatment region for each of the different tissues; constructing an electrode model based on the stimulation electrode array and the constructed target tissue model; and discretizing the target tissue model and the electrode model into a large number of voxels to form a finite element model, wherein the lead field matrix is calculated based on the finite element model.
 14. The apparatus of claim 12, wherein the high-resolution structural images comprise MRI or CT images.
 15. The apparatus of claim 12, wherein the target treatment region comprises a spinal cord or spinal nerves, or a region of the brain.
 16. The apparatus of claim 12, wherein the target treatment region comprises an internal organ.
 17. The apparatus of claim 12, wherein the determining a safety limit for transcutaneous stimulation of the target treatment region by the stimulation electrode array further comprises determining a ratio between current to be applied to the target treatment region and to an avoidance region.
 18. The apparatus of claim 12, wherein the instructions further perform one or more steps comprising: utilizing a stimulation optimization model to provide a feasible solution that includes both the high intensity and focal accuracy of electrical stimulation within the safety limit, wherein the stimulation optimization model is based on the focality of stimulation defining a ratio between the total number of voxels and the number of targeted voxels and the directional intensity of stimulation in the desired direction, wherein the focality and directional intensity of stimulation are balanced by assigning the relative weight terms.
 19. The apparatus of claim 12, wherein the instructions further perform steps comprising: measuring one or more physiological signals using one or more sensors; adjusting one or more changes to the set of stimulation parameters in response to the measured physiological signals; and generating an updated set of stimulation parameters.
 20. The apparatus of claim 19, wherein said one or more physiological signals are selected from a group of physiological signals consisting of EEG, EMG, EKG, ECoG, LFP, spike recording, accelerometer, PPG, SpO2, respiration, and skin impedance signals.
 21. The apparatus of claim 12, wherein the stimulation has a set of stimulation parameters that provide direct current stimulation or alternative current stimulation at a short period of time for focalized stimulation of the target tissue region at any orientation with high precision.
 22. The apparatus of claim 12, wherein the instructions further perform steps comprising designating, simultaneously, at least one target within the target treatment region as a stimulation target and at least one region within the target treatment region as a region to be avoided.
 23. The apparatus of claim 12: wherein the stimulation electrode array comprises multiple stimulation electrodes; and wherein the external device provides independent stimulation signals to each of the multiple electrodes in the stimulation electrode array at a plurality of frequencies, intensities, and/or waveforms.
 24. The apparatus of claim 12, further comprising: one or more sensors selected from a group of physiological sensor consisting of sensors for EEG, EMG, EKG, ECoG, LFP, spike recording, accelerometer, PPG, SpO2, respiration, and skin impedance signals; wherein the external device uses the one or more sensors for monitoring patient conditions, as well as for feedback control and/or feedforward control.
 25. The apparatus of claim 23, wherein the independent stimulation signals are configured to evoke a 2-D electric field pattern in-vitro by selectively choosing the stimulation channels, intensities, polarities, and return channel of selected individual electrodes within the stimulation electrode array.
 26. An apparatus comprising: a stimulation electrode array configured to be placed within a target treatment region of a patient; and an external device comprising: a processor; and a non-transitory memory storing instructions executable by the processor to perform steps comprising: receive high-resolution structural images of the target treatment region; generate 3D model based on the high-resolution structural images of the target treatment region with the stimulation electrode array placed within the target treatment region; calculate, based on the 3D model, a lead field matrix associated with the stimulation electrode array and the target treatment region; determine a safety limit for transcutaneous stimulation of the target treatment region by the stimulation electrode array; and generate a set of stimulation parameters to be delivered by the stimulation electrode array based on the lead field matrix and the safety limit by: segmenting the high-resolution structural images into different tissues according to grey levels; generating a target tissue model of the target treatment region for each of the different tissues; constructing an electrode model based on the stimulation electrode array and the constructed target tissue model; and discretizing the target tissue model and the electrode model into a large number of voxels to form a finite element model, wherein the lead field matrix is calculated based on the finite element model; and communicate to the stimulation electrode array the set of stimulation parameters; wherein the stimulation electrode array is configured to deliver a stimulation with high intensity and focal accuracy of electrical stimulation within the safety limit based on the set of stimulation parameters, where the set of stimulation parameters assigns a relative weight to a directional intensity and a focality of stimulation of the target treatment region. 